*Figure 1 and A14
foreach var in std_missnojd std_olust wpch_ee_d3 leave_cons {
	foreach sample in  man_main woman_main{
	preserve
	local s=0
	matrix C =J(1,3,.)
	local s=`s'+1
	reg `var' hara_comb  [pw=gen_edu2_alder5_w ] if `sample' ==1
	matrix C[1,1]=_b[ hara_comb]
	matrix C[1,2]=(_b[ hara_comb]+1.96*_se[ hara_comb])
	matrix C[1,3]=(_b[ hara_comb]-1.96*_se[ hara_comb])

	svmat C
	ren C1 est_`s'
	ren C2 ci_low_`s'
	ren C3 ci_high_`s'
	gen pos_`s'=8 if est_`s'!=.
	
	*****
	matrix C =J(1,3,.)
	local s=`s'+1
	xi: areg `var' hara_comb   [pw=gen_edu2_alder5_w ] if `sample' ==1 ,  abs(ssyk3)
	matrix C[1,1]=_b[ hara_comb]
	matrix C[1,2]=(_b[ hara_comb]+1.96*_se[ hara_comb])
	matrix C[1,3]=(_b[ hara_comb]-1.96*_se[ hara_comb])

	svmat C
	ren C1 est_`s'
	ren C2 ci_low_`s'
	ren C3 ci_high_`s'
	gen pos_`s'=7 if est_`s'!=.
	
	****
	local s=`s'+1
	xi: reg `var' hara_comb  fem_bos  [pw=gen_edu2_alder5_w ] if `sample' ==1 ,
	matrix C[1,1]=_b[ hara_comb]
	matrix C[1,2]=(_b[ hara_comb]+1.96*_se[ hara_comb])
	matrix C[1,3]=(_b[ hara_comb]-1.96*_se[ hara_comb])

	svmat C
	ren C1 est_`s'
	ren C2 ci_low_`s'
	ren C3 ci_high_`s'
	gen pos_`s'=6 if est_`s'!=.
	
	local s=`s'+1
	xi: reg `var' hara_comb i.age_cat i.edu3 i.famtyp europe non_europe  [pw=gen_edu2_alder5_w ] if `sample' ==1 , 
	matrix C[1,1]=_b[ hara_comb]
	matrix C[1,2]=(_b[ hara_comb]+1.96*_se[ hara_comb])
	matrix C[1,3]=(_b[ hara_comb]-1.96*_se[ hara_comb])

	svmat C
	ren C1 est_`s'
	ren C2 ci_low_`s'
	ren C3 ci_high_`s'
	gen pos_`s'=5 if est_`s'!=.
	****
	local s=`s'+1
	xi: reg `var' hara_comb  workpl_m [pw=gen_edu2_alder5_w ] if `sample' ==1 , 
	matrix C[1,1]=_b[ hara_comb]
	matrix C[1,2]=(_b[ hara_comb]+1.96*_se[ hara_comb])
	matrix C[1,3]=(_b[ hara_comb]-1.96*_se[ hara_comb])

	svmat C
	ren C1 est_`s'
	ren C2 ci_low_`s'
	ren C3 ci_high_`s'
	gen pos_`s'=4 if est_`s'!=.
	****
	local s=`s'+1
	xi: reg `var' hara_comb pc1 [pw=gen_edu2_alder5_w ] if `sample' ==1 ,   
	matrix C[1,1]=_b[ hara_comb]
	matrix C[1,2]=(_b[ hara_comb]+1.96*_se[ hara_comb])
	matrix C[1,3]=(_b[ hara_comb]-1.96*_se[ hara_comb])

	svmat C
	ren C1 est_`s'
	ren C2 ci_low_`s'
	ren C3 ci_high_`s'
	gen pos_`s'=3 if est_`s'!=.

	****
	local s=`s'+1
	xi: areg `var' hara_comb pc1  i.age_cat fem_bos i.edu3 i.famtyp europe non_europe workpl_m if `sample' ==1  [pw=gen_edu2_alder5_w ] , abs(ssyk3)
	matrix C[1,1]=_b[ hara_comb]
	matrix C[1,2]=(_b[ hara_comb]+1.96*_se[ hara_comb])
	matrix C[1,3]=(_b[ hara_comb]-1.96*_se[ hara_comb])
	svmat C
	ren C1 est_`s'
	ren C2 ci_low_`s'
	ren C3 ci_high_`s'
	gen pos_`s'=2 if est_`s'!=.
	
	local s=`s'+1
	xi: areg `var' hara_comb if `sample' ==1  [pw=gen_edu2_alder5_w ] , abs(cfarnr_lopnr)
	matrix C[1,1]=_b[ hara_comb]
	matrix C[1,2]=(_b[ hara_comb]+1.96*_se[ hara_comb])
	matrix C[1,3]=(_b[ hara_comb]-1.96*_se[ hara_comb])

	svmat C
	ren C1 est_`s'
	ren C2 ci_low_`s'
	ren C3 ci_high_`s'
	gen pos_`s'=1 if est_`s'!=.



	* top part-yscale
	twoway (scatter  pos_1 est_1  , msize(medlarge)  msymbol(circle)   mcolor(gs1)) ///
	(rbar ci_low_1 ci_high_1 pos_1 , lwidth(vvthin)  fcolor(gs1) lcolor(gs1) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_2 est_2  , msize(medlarge)  msymbol(circle)  mcolor(gs3)) ///
	(rbar ci_low_2 ci_high_2 pos_2 , lwidth(vvthin)  fcolor(gs3) lcolor(gs3) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_3 est_3  , msize(medlarge)  msymbol(circle)  mcolor(gs3)) ///
	(rbar ci_low_3 ci_high_3 pos_3 , lwidth(vvthin)  fcolor(gs5) lcolor(gs5) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_4 est_4  , msize(medlarge)  msymbol(circle)  mcolor(gs5)) ///
	(rbar ci_low_4 ci_high_4 pos_4 , lwidth(vvthin)  fcolor(gs7) lcolor(gs7) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_5 est_5  , msize(medlarge)  msymbol(circle)  mcolor(gs9)) ///
	(rbar ci_low_5 ci_high_5 pos_5 , lwidth(vvthin)  fcolor(gs9) lcolor(gs9) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_6 est_6  , msize(medlarge)  msymbol(circle)  mcolor(gs11)) ///
	(rbar ci_low_6 ci_high_6 pos_6 , lwidth(vvthin)  fcolor(gs11) lcolor(gs11) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_7 est_7  , msize(medlarge)  msymbol(circle)  mcolor(gs13)) ///
	(rbar ci_low_7 ci_high_7 pos_7 , lwidth(vvthin)  fcolor(gs13) lcolor(gs13) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_8 est_8  , msize(medlarge)  msymbol(square)  mcolor(gs1)) ///
	(rbar ci_low_8 ci_high_8 pos_8 , lwidth(vvthin)  fcolor(gs1) lcolor(gs1) lpattern(solid) horizontal barwidth(.015)) ///
	,xline(0, lcolor(gs7) lpattern(dash)) yscale(range(0 9)) xlabel(,grid) ///
	 ytitle(, size(zero)) ylabel(none)  scale(.7) ysize(20) xsize(20) scheme(s1mono) title(`sample' `var') ///
	 legend(on order(1 "No controls" 3 "Occupation Fixed Effects" 5 "Female Supervisor" 7 "Demographics"  ///
	 9 "Workplace Share of Men"  11 "PCA Bad Work Conditions"  13  "All Controls" 15 "Occupation Fixed Effects" )  col(3)) 
	graph save reg_`sample'_`var', replace
	restore
	
	}
}
		
				
graph combine reg_woman_main_std_missnojd.gph reg_woman_main_std_olust.gph ///
	reg_man_main_std_missnojd.gph  reg_man_main_std_olust.gph , ///
	ysize(18) xsize(20) col(2) row(2) iscale(1) scheme(s1mono) title(Women Men)
graph save Figure_1, replace
		
graph combine reg_woman_main_wpch_ee_d3.gph reg_man_main_wpch_ee_d3.gph ///
	reg_woman_main_leave_cons.gph reg_man_main_leave_cons.gph , ///
	ysize(18) xsize(20) col(2) row(2) iscale(1) scheme(s1mono) title(Women Men)
graph save Figure A14, replace
		

/****************************************************/
*define program to make graph with binned averages in a specific number of bins*/
/****************************************************/

capture program drop bin_split_line
program bin_split_line

	*******the arguments are the outcome variable, the x variable, weight, restriction on the sample used, first subsample, second subsample, Y axis title, X axis, name on first subsample, name second subsample, graph title**************
	args ss vs  bin weight lim lim1 lim2 tit1 tit2  name1 name2 grtitle
		
	preserve
	********define ranks to construct bins****************
	xtile bin1= `vs' [pw=`weight'] if `lim' &`ss'!=. &`vs'!=. & `lim1', nq(`bin')
	xtile bin2= `vs' [pw=`weight'] if `lim' &`ss'!=. &`vs'!=. & `lim2', nq(`bin')
	************create bin averages *************
	bysort bin1: egen mean`vs'1=mean(`vs') if  bin1!=. & `lim' & `lim1'
	bysort bin1: egen mean`ss'1=mean(`ss') if `ss'!=. & bin1!=. & `lim' & `lim1'
	bysort bin2: egen mean`vs'2=mean(`vs') if  bin2!=. & `lim' & `lim2'
	bysort bin2: egen mean`ss'2=mean(`ss') if `ss'!=. & bin2!=. & `lim' & `lim2'
	egen tag1=tag(bin1) if  bin1!=. & `lim' & `lim1'
	egen tag2=tag(bin2) if  bin2!=. & `lim' & `lim2'
	
	egen max1 =max(mean`vs'1  )
	egen min1 =min(mean`vs'1 )
    egen max2 =max( mean`vs'2 )
	egen min2 =min( mean`vs'2)
	egen max=rowmax(max1 max2)
	egen min=rowmin(min1 min2)
	twoway  (scatter mean`ss'1 mean`vs'1 if tag1==1 ,  msize(medium) msymbol(circle_hollow) mcolor(gs6)) ///
	(scatter mean`ss'2 mean`vs'2 if tag2==1 , msize(medium) msymbol(triangle) mcolor(black)) ///
	(qfit `ss' `vs' if  `lim' & `lim1' &`vs'>min & `vs'<max ,  lcolor(gs6)) ///
	(qfit `ss' `vs' if  `lim' & `lim2'  &`vs'>min & `vs'<max  , lcolor(black)), ///
    title(`grtitle')   ytitle(`tit1') xtitle(`tit2') 	ylabel(#5, angle(horizontal) grid)   xlabel(#5) legend(on order(1 "`name1'" 2 "`name2'")) scheme(s1mono)
restore

end

*Right hand side of Figure 2, A4 and A6
foreach var in   hara_comb hara_coll sex_trak high_freq vhigh_freq{
	foreach treat in workpl_m{
	preserve
	foreach sample in  man woman{

	local s=0
	matrix C =J(1,3,.)
	local s=`s'+1
	reg `var' `treat'  [pw=gen_edu2_alder5_w ] if main_samp==1 & `sample' ==1
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=7 if est_`s'_`sample'!=.
	
	*****
	matrix C =J(1,3,.)
	local s=`s'+1
	xi: areg `var' `treat' [pw=gen_edu2_alder5_w ] if  main_samp==1 & `sample' ==1  , abs(ssyk3) 
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=6 if est_`s'_`sample'!=.
	
	****
	local s=`s'+1
	xi: reg `var' `treat' fem_bos  [pw=gen_edu2_alder5_w ] if  main_samp==1 & `sample' ==1 ,
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=5 if est_`s'_`sample'!=.
	****
	local s=`s'+1
	xi: reg `var' `treat' i.age_cat i.edu3 i.famtyp europe non_europe  [pw=gen_edu2_alder5_w ] if  main_samp==1 & `sample' ==1
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=4 if est_`s'_`sample'!=.
	****
	local s=`s'+1
	reg `var' `treat' high_ind_fe [pw=gen_edu2_alder5_w ] if  main_samp==1 & `sample' ==1 , 
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=3 if est_`s'_`sample'!=.

	local s=`s'+1
	areg `var' `treat'[pw=gen_edu2_alder5_w ] if  main_samp==1 & `sample' ==1 ,  abs( fsize_cat) 
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=2 if est_`s'_`sample'!=.
	local s=`s'+1
	xi: reghdfe `var' `treat' i.fsize_cat high_ind_fe i.age_cat fem_bos i.edu3 i.famtyp europe non_europe [pw=gen_edu2_alder5_w ] if  main_samp==1 & `sample' ==1 , abs(ssyk3)
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=1 if est_`s'_`sample'!=.
}

	twoway (scatter  pos_1_w est_1_w  , msize(medlarge)  msymbol(circle_hollow)   mcolor(gs1)) ///
	(rbar ci_low_1_w ci_high_1_w pos_1_w , lwidth(vvthin)  fcolor(gs1) lcolor(gs1) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_2_w est_2_w  , msize(medlarge)  msymbol(circle_hollow)  mcolor(gs4)) ///
	(rbar ci_low_2_w ci_high_2_w pos_2_w , lwidth(vvthin)  fcolor(gs4) lcolor(gs4) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_3_w est_3_w  , msize(medlarge)  msymbol(circle_hollow)  mcolor(gs6)) ///
	(rbar ci_low_3_w ci_high_3_w pos_3_w , lwidth(vvthin)  fcolor(gs6) lcolor(gs6) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_4_w est_4_w  , msize(medlarge)  msymbol(circle_hollow)  mcolor(gs8)) ///
	(rbar ci_low_4_w ci_high_4_w pos_4_w , lwidth(vvthin)  fcolor(gs8) lcolor(gs8) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_5_w est_5_w  , msize(medlarge)  msymbol(circle_hollow)  mcolor(gs10)) ///
	(rbar ci_low_5_w ci_high_5_w pos_5_w , lwidth(vvthin)  fcolor(gs10) lcolor(gs10) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_6_w est_6_w  , msize(medlarge)  msymbol(circle_hollow)  mcolor(gs12)) ///
	(rbar ci_low_6_w ci_high_6_w pos_6_w , lwidth(vvthin)  fcolor(gs12) lcolor(gs12) lpattern(solid) horizontal barwidth(.015)) ///	
	(scatter  pos_7_w est_7_w  , msize(medlarge)  msymbol(circle_hollow)  mcolor(gs14)) ///
	(rbar ci_low_7_w ci_high_7_w pos_7_w , lwidth(vvthin)  fcolor(gs14) lcolor(gs14) lpattern(solid) horizontal barwidth(.015)) ///	
	(scatter  pos_1_m est_1_m  , msize(medlarge)  msymbol(triangle)   mcolor(gs1)) ///
	(rbar ci_low_1_m ci_high_1_m pos_1_m , lwidth(vvthin)  fcolor(gs1) lcolor(gs1) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_2_m est_2_m  , msize(medlarge)  msymbol(triangle)  mcolor(gs4)) ///
	(rbar ci_low_2_m ci_high_2_m pos_2_m , lwidth(vvthin)  fcolor(gs3) lcolor(gs4) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_3_m est_3_m  , msize(medlarge)  msymbol(triangle)  mcolor(gs6)) ///
	(rbar ci_low_3_m ci_high_3_m pos_3_m , lwidth(vvthin)  fcolor(gs6) lcolor(gs6) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_4_m est_4_m  , msize(medlarge)  msymbol(triangle)  mcolor(gs6)) ///
	(rbar ci_low_4_m ci_high_4_m pos_4_m , lwidth(vvthin)  fcolor(gs8) lcolor(gs8) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_5_m est_5_m  , msize(medlarge)  msymbol(triangle)  mcolor(gs10)) ///
	(rbar ci_low_5_m ci_high_5_m pos_5_m , lwidth(vvthin)  fcolor(gs10) lcolor(gs10) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_6_m est_6_m  , msize(medlarge)  msymbol(triangle)  mcolor(gs12)) ///
	(rbar ci_low_6_m ci_high_6_m pos_6_m , lwidth(vvthin)  fcolor(gs12) lcolor(gs12) lpattern(solid) horizontal barwidth(.015)) ///
   	(scatter  pos_7_m est_7_m  , msize(medlarge)  msymbol(triangle)  mcolor(gs14)) ///
	(rbar ci_low_7_m ci_high_7_m pos_7_m , lwidth(vvthin)  fcolor(gs14) lcolor(gs14) lpattern(solid) horizontal barwidth(.015)) ///
	,xline(0, lcolor(gs7) lpattern(dash)) yscale(range(0 8))  xlabel(, grid) ///
	 ytitle(, size(zero)) scale(1) ylabel(none)   ysize(20) xsize(20) scheme(s1mono) title(`sample' `var') ///
	 legend(on order(1 "No controls" 3 "Occupation Fixed Effects" 5 "Female Supervisor" 7 "Demographics"  ///
	 9 "Productivity"  11 "Firm Size"  13  "All Controls" )  col(3)) ///
	 xtitle(Estimated Relationship Betweeen Share of Men in the Workplace and Sexual Harassment)
	graph save reg_`sample'_`treat'_`var', replace

	restore
	}
}

 
 *Left hand side of Figure 2, A4 and A6, and all figures for Figure A5
foreach var in   hara_comb hara_coll sex_trak high_freq vhigh_fre  bully  coll_napp chef_napp  coll_nsup chef_nsup   conf_coll  conf_boss pc1{
	bin_split_line `var' workpl_m  100 gen_edu2_alder5_w  main_sample==1  woman==1 woman==0 "Share Sexually Harassed in the Last 12 Months" ///
		"Share of Men in Workplace"  Women Men
	graph save " `var'_wpm.gph", replace  
}
graph combine hara_wpm.gph reg__workpl_m_hara_comb.gph , col(2) scheme(s1mono)  ysize(7) xsize(20)  iscale(1.4)
graph save Figure_2, replace
	
graph combine sex_trak_wpm.gph reg__workpl_m_sex_trak.gph ///
hara_coll_wpm.gph reg__workpl_m_hara_coll.gph ///
, col(2) scheme(s1mono)  ysize(11) xsize(20)  iscale(1)
graph save Figure_A4, replace	

graph combine bully_wpm.gph conf_coll_wpm.gph conf_boss_wpm.gph ///
	coll_napp_wpm.gph chef_napp_wpm.gph  coll_nsup_wpm.gph chef_nsup_wpm.gph   pc1_wpm.gph, ///
	col(2) scheme(s1mono)  ysize(20) xsize(20)  iscale(.5)
graph save Figure_A5, replace
log using reg_tofigA5, replace
foreach var in  bully conf_coll conf_boss  coll_napp chef_napp  coll_nsup chef_nsup  pc1 {
	reg `var' workpl_m if main==1 & woman==1
	reg `var' workpl_m if main==1 & woman==0
}
		
graph combine high_freq_wpm.gph reg__workpl_m_high_freq.gph ///
vhigh_freq_wpm.gph reg__workpl_m_vhigh_freq.gph ///
, col(2) scheme(s1mono)  ysize(11) xsize(20)  iscale(1)
graph save Figure_A6, replace
				
**Figure A3
bin_split_line hara_comb workpl_m  50 gen_edu2_alder5_w  "main_sample==1 & cfarnr_lopnr!=. & woman==1"  edu3<3 edu3==3 "Share Sexually Harassed in the Last 12 Months" ///
"Share of Men in the Workplace" "Below Tertiary Education" "Tertiary Education"
graph save hara_wpm_hle_w.gph, replace
bin_split_line hara_comb workpl_m  50 gen_edu2_alder5_w  "main_sample==1 & cfarnr_lopnr!=. & woman==0"  edu3<3 edu3==3 "Share Sexually Harassed in the Last 12 Months" ///
"Share of Men in the Workplace"  "Below Tertiary Education" "Tertiary Education"
graph save hara_wpm_hle_m.gph, replace
graph combine hara_wpm_hle_w.gph hara_wpm_hle_m.gph , col(2) scheme(s1mono)  ysize(7) xsize(20)  iscale(1.4)
graph save Figure_A3, replace		
	
	
*Figure A13			
bin_split_line hara_comb workpl_m  50 gen_edu2_alder5_w  "main_sample==1 & cfarnr_lopnr!=. & woman==1"  hw==1 hw==0 "Share Sexually Harassed in the Last 12 Months" ///
"Share of Men in the Workplace"  "High Wage Workplace" "Low Wage Workplace"
graph save hara_wpm_hlw_w.gph, replace			
bin_split_line hara_comb workpl_m  50 gen_edu2_alder5_w  "main_sample==1 & cfarnr_lopnr!=. & woman==0"  hw==1 hw==0 "Share Sexually Harassed in the Last 12 Months" ///
"Share of Men in the Workplace"  "High Wage Workplace" "Low Wage Workplace"
raph save hara_wpm_hlw_m.gph, replace

graph combine hara_wpm_hlw_w.gph hara_wpm_hlw_m.gph , col(2) scheme(s1mono)  ysize(7) xsize(20)  iscale(1.4)
graph save Figure_A13, replace
	
***Figure A2
foreach var in   hara_comb  {
foreach treat in ssyk3_m{
	preserve
	foreach sample in  man woman{

	local s=0
	matrix C =J(1,3,.)
	local s=`s'+1
	reg `var' `treat'  [pw=gen_edu2_alder5_w ] if main_samp==1 & `sample' ==1
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=7 if est_`s'_`sample'!=.
	
	*****
	matrix C =J(1,3,.)
	local s=`s'+1
	xi: reg `var' `treat' workpl_m [pw=gen_edu2_alder5_w ] if  main_samp==1 & `sample' ==1  , 
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=6 if est_`s'_`sample'!=.
	
	****
	local s=`s'+1
	xi: reg `var' `treat' fem_bos  [pw=gen_edu2_alder5_w ] if  main_samp==1 & `sample' ==1 ,
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=5 if est_`s'_`sample'!=.
	****
	local s=`s'+1
	xi: reg `var' `treat' i.age_cat i.edu3 i.famtyp europe non_europe  [pw=gen_edu2_alder5_w ] if  main_samp==1 & `sample' ==1
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=4 if est_`s'_`sample'!=.
	****
	local s=`s'+1
	reg `var' `treat' high_ind_fe [pw=gen_edu2_alder5_w ] if  main_samp==1 & `sample' ==1 , 
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=3 if est_`s'_`sample'!=.

	local s=`s'+1
	areg `var' `treat'[pw=gen_edu2_alder5_w ] if  main_samp==1 & `sample' ==1 ,  abs( fsize_cat) 
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=2 if est_`s'_`sample'!=.
	local s=`s'+1
	xi: reg `var' `treat' i.fsize_cat  i.age_cat fem_bos i.edu3 i.famtyp workpl_m europe non_europe [pw=gen_edu2_alder5_w ] if  main_samp==1 & `sample' ==1 , 
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=1 if est_`s'_`sample'!=.

	}

twoway (scatter  pos_1_w est_1_w  , msize(medlarge)  msymbol(circle_hollow)   mcolor(gs1)) ///
	(rbar ci_low_1_w ci_high_1_w pos_1_w , lwidth(vvthin)  fcolor(gs1) lcolor(gs1) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_2_w est_2_w  , msize(medlarge)  msymbol(circle_hollow)  mcolor(gs4)) ///
	(rbar ci_low_2_w ci_high_2_w pos_2_w , lwidth(vvthin)  fcolor(gs4) lcolor(gs4) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_3_w est_3_w  , msize(medlarge)  msymbol(circle_hollow)  mcolor(gs6)) ///
	(rbar ci_low_3_w ci_high_3_w pos_3_w , lwidth(vvthin)  fcolor(gs6) lcolor(gs6) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_4_w est_4_w  , msize(medlarge)  msymbol(circle_hollow)  mcolor(gs8)) ///
	(rbar ci_low_4_w ci_high_4_w pos_4_w , lwidth(vvthin)  fcolor(gs8) lcolor(gs8) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_5_w est_5_w  , msize(medlarge)  msymbol(circle_hollow)  mcolor(gs10)) ///
	(rbar ci_low_5_w ci_high_5_w pos_5_w , lwidth(vvthin)  fcolor(gs10) lcolor(gs10) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_6_w est_6_w  , msize(medlarge)  msymbol(circle_hollow)  mcolor(gs12)) ///
	(rbar ci_low_6_w ci_high_6_w pos_6_w , lwidth(vvthin)  fcolor(gs12) lcolor(gs12) lpattern(solid) horizontal barwidth(.015)) ///	
	(scatter  pos_7_w est_7_w  , msize(medlarge)  msymbol(circle_hollow)  mcolor(gs14)) ///
	(rbar ci_low_7_w ci_high_7_w pos_7_w , lwidth(vvthin)  fcolor(gs14) lcolor(gs14) lpattern(solid) horizontal barwidth(.015)) ///	
	(scatter  pos_1_m est_1_m  , msize(medlarge)  msymbol(triangle)   mcolor(gs1)) ///
	(rbar ci_low_1_m ci_high_1_m pos_1_m , lwidth(vvthin)  fcolor(gs1) lcolor(gs1) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_2_m est_2_m  , msize(medlarge)  msymbol(triangle)  mcolor(gs4)) ///
	(rbar ci_low_2_m ci_high_2_m pos_2_m , lwidth(vvthin)  fcolor(gs3) lcolor(gs4) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_3_m est_3_m  , msize(medlarge)  msymbol(triangle)  mcolor(gs6)) ///
	(rbar ci_low_3_m ci_high_3_m pos_3_m , lwidth(vvthin)  fcolor(gs6) lcolor(gs6) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_4_m est_4_m  , msize(medlarge)  msymbol(triangle)  mcolor(gs6)) ///
	(rbar ci_low_4_m ci_high_4_m pos_4_m , lwidth(vvthin)  fcolor(gs8) lcolor(gs8) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_5_m est_5_m  , msize(medlarge)  msymbol(triangle)  mcolor(gs10)) ///
	(rbar ci_low_5_m ci_high_5_m pos_5_m , lwidth(vvthin)  fcolor(gs10) lcolor(gs10) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_6_m est_6_m  , msize(medlarge)  msymbol(triangle)  mcolor(gs12)) ///
	(rbar ci_low_6_m ci_high_6_m pos_6_m , lwidth(vvthin)  fcolor(gs12) lcolor(gs12) lpattern(solid) horizontal barwidth(.015)) ///
   	(scatter  pos_7_m est_7_m  , msize(medlarge)  msymbol(triangle)  mcolor(gs14)) ///
	(rbar ci_low_7_m ci_high_7_m pos_7_m , lwidth(vvthin)  fcolor(gs14) lcolor(gs14) lpattern(solid) horizontal barwidth(.015)) ///
   ,xline(0, lcolor(gs7) lpattern(dash)) yscale(range(0 8))  xlabel(, grid) ///
	 ytitle(, size(zero)) scale(1) ylabel(none)   ysize(20) xsize(20) scheme(s1mono) title(`sample' `var') ///
	 legend(on order(1 "No controls" 3 "WP share of men" 5 "Female Supervisor" 7 "Demographics"  ///
	 9 "Productivity"  11 "Firm Size"  13  "All Controls" )  col(3)) ///
	 xtitle(Estimated Relationship Betweeen Share of Men in the Occupation and Sexual Harassment)
	graph save reg_`sample'_`treat'_`var', replace
	restore
	}	
}
bin_split_line hara_comb ssyk3_m 100 gen_edu2_alder5_w  main_sample==1  woman==1 woman==0 "Share Sexually Harassed in the Last 12 Months" ///
graph save "hara_comb_occm.gph", replace  
graph save Figure_A2, replace

*Figure A7 alternative weights
bin_split_line `var' workpl_m  100 gen_edu2_alder3_pinc3_br_w  main_sample==1  woman==1 woman==0 "Share Sexually Harassed in the Last 12 Months" ///
"Share of Men in Workplace"  Women Men
graph save "hara_comb_wpm_aw.gph", replace  

*Right hand side of Figure A7 with alternative weights
foreach var in   hara_comb {
	foreach treat in workpl_m{
	preserve
	foreach sample in  man woman{

	local s=0
	matrix C =J(1,3,.)
	local s=`s'+1
	reg `var' `treat'  [pw=gen_edu2_alder3_pinc3_br_w ] if main_samp==1 & `sample' ==1
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=7 if est_`s'_`sample'!=.
	
	*****
	matrix C =J(1,3,.)
	local s=`s'+1
	xi: areg `var' `treat' [pw=gen_edu2_alder3_pinc3_br_w ] if  main_samp==1 & `sample' ==1  , abs(ssyk3) 
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=6 if est_`s'_`sample'!=.
	
	****
	local s=`s'+1
	xi: reg `var' `treat' fem_bos  [pw=gen_edu2_alder3_pinc3_br_w ] if  main_samp==1 & `sample' ==1 ,
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=5 if est_`s'_`sample'!=.
	****
	local s=`s'+1
	xi: reg `var' `treat' i.age_cat i.edu3 i.famtyp europe non_europe  [pw=gen_edu2_alder3_pinc3_br_w ] if  main_samp==1 & `sample' ==1
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=4 if est_`s'_`sample'!=.
	****
	local s=`s'+1
	reg `var' `treat' high_ind_fe [pw=gen_edu2_alder3_pinc3_br_w ] if  main_samp==1 & `sample' ==1 , 
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=3 if est_`s'_`sample'!=.

	local s=`s'+1
	areg `var' `treat'[pw=gen_edu2_alder3_pinc3_br_w ] if  main_samp==1 & `sample' ==1 ,  abs( fsize_cat) 
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=2 if est_`s'_`sample'!=.
	local s=`s'+1
	xi: reghdfe `var' `treat' i.fsize_cat high_ind_fe i.age_cat fem_bos i.edu3 i.famtyp europe non_europe [pw=gen_edu2_alder3_pinc3_br_w ] if  main_samp==1 & `sample' ==1 , abs(ssyk3)
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'_`sample'
	ren C2 ci_low_`s'_`sample'
	ren C3 ci_high_`s'_`sample'
	gen pos_`s'_`sample'=1 if est_`s'_`sample'!=.

}

twoway (scatter  pos_1_w est_1_w  , msize(medlarge)  msymbol(circle_hollow)   mcolor(gs1)) ///
	(rbar ci_low_1_w ci_high_1_w pos_1_w , lwidth(vvthin)  fcolor(gs1) lcolor(gs1) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_2_w est_2_w  , msize(medlarge)  msymbol(circle_hollow)  mcolor(gs4)) ///
	(rbar ci_low_2_w ci_high_2_w pos_2_w , lwidth(vvthin)  fcolor(gs4) lcolor(gs4) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_3_w est_3_w  , msize(medlarge)  msymbol(circle_hollow)  mcolor(gs6)) ///
	(rbar ci_low_3_w ci_high_3_w pos_3_w , lwidth(vvthin)  fcolor(gs6) lcolor(gs6) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_4_w est_4_w  , msize(medlarge)  msymbol(circle_hollow)  mcolor(gs8)) ///
	(rbar ci_low_4_w ci_high_4_w pos_4_w , lwidth(vvthin)  fcolor(gs8) lcolor(gs8) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_5_w est_5_w  , msize(medlarge)  msymbol(circle_hollow)  mcolor(gs10)) ///
	(rbar ci_low_5_w ci_high_5_w pos_5_w , lwidth(vvthin)  fcolor(gs10) lcolor(gs10) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_6_w est_6_w  , msize(medlarge)  msymbol(circle_hollow)  mcolor(gs12)) ///
	(rbar ci_low_6_w ci_high_6_w pos_6_w , lwidth(vvthin)  fcolor(gs12) lcolor(gs12) lpattern(solid) horizontal barwidth(.015)) ///	
	(scatter  pos_7_w est_7_w  , msize(medlarge)  msymbol(circle_hollow)  mcolor(gs14)) ///
	(rbar ci_low_7_w ci_high_7_w pos_7_w , lwidth(vvthin)  fcolor(gs14) lcolor(gs14) lpattern(solid) horizontal barwidth(.015)) ///	
	(scatter  pos_1_m est_1_m  , msize(medlarge)  msymbol(triangle)   mcolor(gs1)) ///
	(rbar ci_low_1_m ci_high_1_m pos_1_m , lwidth(vvthin)  fcolor(gs1) lcolor(gs1) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_2_m est_2_m  , msize(medlarge)  msymbol(triangle)  mcolor(gs4)) ///
	(rbar ci_low_2_m ci_high_2_m pos_2_m , lwidth(vvthin)  fcolor(gs3) lcolor(gs4) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_3_m est_3_m  , msize(medlarge)  msymbol(triangle)  mcolor(gs6)) ///
	(rbar ci_low_3_m ci_high_3_m pos_3_m , lwidth(vvthin)  fcolor(gs6) lcolor(gs6) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_4_m est_4_m  , msize(medlarge)  msymbol(triangle)  mcolor(gs6)) ///
	(rbar ci_low_4_m ci_high_4_m pos_4_m , lwidth(vvthin)  fcolor(gs8) lcolor(gs8) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_5_m est_5_m  , msize(medlarge)  msymbol(triangle)  mcolor(gs10)) ///
	(rbar ci_low_5_m ci_high_5_m pos_5_m , lwidth(vvthin)  fcolor(gs10) lcolor(gs10) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_6_m est_6_m  , msize(medlarge)  msymbol(triangle)  mcolor(gs12)) ///
	(rbar ci_low_6_m ci_high_6_m pos_6_m , lwidth(vvthin)  fcolor(gs12) lcolor(gs12) lpattern(solid) horizontal barwidth(.015)) ///
   	(scatter  pos_7_m est_7_m  , msize(medlarge)  msymbol(triangle)  mcolor(gs14)) ///
	(rbar ci_low_7_m ci_high_7_m pos_7_m , lwidth(vvthin)  fcolor(gs14) lcolor(gs14) lpattern(solid) horizontal barwidth(.015)) ///
   ,xline(0, lcolor(gs7) lpattern(dash)) yscale(range(0 8))  xlabel(, grid) ///
	 ytitle(, size(zero)) scale(1) ylabel(none)   ysize(20) xsize(20) scheme(s1mono) title(`sample' `var') ///
	 legend(on order(1 "No controls" 3 "Occupation Fixed Effects" 5 "Female Supervisor" 7 "Demographics"  ///
	 9 "Productivity"  11 "Firm Size"  13  "All Controls" )  col(3)) ///
	 xtitle(Estimated Relationship Betweeen Share of Men in the Workplace and Sexual Harassment)
	graph save reg_`sample'_`treat'_`var'_aw, replace
	restore
	}
}
graph combine hara_wpm_aw.gph reg__workpl_m_hara_comb_aw.gph , col(2) scheme(s1mono)  ysize(7) xsize(20)  iscale(1.4)
graph save Figure A7, replace

*Figure V	
* define program to make graph with seat share and vote share zero bin and fixed observations in bins*/
/****************************************************/
capture program drop haragaput_wage_shm
program haragaput_wage_shm

	*******the arguments are the outcome variable, the forcing variable, restriction on the sample used, ***
	*****restriction on the distance from the threshold, Y axis title, X axis title**************
	args ss ss2 vs  nbin weight lim lim1 lim2  ut name1 name2 
	sum `vs'
	preserve
	********define bins****************
	gen bin= int(`vs'*`nbin') / `nbin'
	gen xvar=bin+.5/ `nbin'
	************create bin averages *************
	bysort bin: egen y1_s1b=mean(`ss') if `ss'!=. & bin!=. & `lim' & `lim1'
	bysort bin: egen y1_s1=max(y1_s1b) if `ss'!=. & bin!=. & `lim' 
	bysort bin: egen y1_s2b=mean(`ss') if `ss'!=. & bin!=. & `lim' & `lim2'
	bysort bin: egen y1_s2=max(y1_s2b) if `ss'!=. & bin!=. & `lim' 
	gen y1_diff =(y1_s1-y1_s2)*`ut'
    bysort bin: egen y2=mean(`ss2') if `ss2'!=. & bin!=. & `lim' 
	egen tag_y1_diff=tag(bin) if  bin!=. & y1_diff!=.
	egen tag_y2=tag(bin) if  bin!=. & y2!=.
	twoway  (scatter  y1_diff xvar if tag_y1_diff==1 ,  yaxis(1) msize(large) msymbol(circle_hollow) mcolor(black)) ///
    (scatter  y2 xvar if tag_y2==1 ,  yaxis(2) msize(medium) msymbol(square) mcolor(gs6)), ///
	xlabel(#5) legend(on order(1 "`name1'" 2 "`name2'" )) scheme(s1mono) ///
  	ylabel(#5, angle(horizontal) grid)
end 

haragaput_wage_shm hara_comb firm_fe_5 workpl_m 20 gen_edu2_alder5_w ///
	"main_samp==1 & firm_fe_5!=." "woman==1" "woman==0""-.11" 
	
graph save  Figure_5, replace
haragaput_wage_shm hara_comb wp_wage_level workpl_m 20 gen_edu2_alder5_w ///
	"main_samp==1 & firm_fe_5!=." "woman==1" "woman==0""-.11" 
graph save  Figure_A11, replace

*Generate Figure A12
/****************************************************/
capture program drop var3_shm
program var3_shm

	*******the arguments are the outcome variables, the x variable, restriction on the sample used, ***
	*****restriction on the distance from the threshold, Y axis title, X axis title**************
	args ss ss2 ss3 vs  nbin lim  name1 name2   name3
	preserve
	replace `vs' =`vs'-.0001 if `vs'==1
	********define rabins****************
	gen bin= int(`vs'*`nbin') / `nbin'
	gen xvar=bin+.5/ `nbin'
	************create bin averages of both forcing variable and outcome*************
	bysort bin: egen y1=mean(`ss') if `ss'!=. & bin!=. & `lim' 
    bysort bin: egen y2=mean(`ss2') if `ss2'!=. & bin!=. & `lim' 
	bysort bin: egen y3=mean(`ss3') if `ss3'!=. & bin!=. & `lim' 
	egen  tag_y1=tag(bin) if  bin!=. & y1!=.
	egen tag_y2=tag(bin) if  bin!=. & y2!=.
    egen tag_y3=tag(bin) if  bin!=. & y3!=.
	twoway (line y1 xvar if tag_y1==1 ,sort lpattern(solid) lcolor(black)) ///
    (line y2 xvar if tag_y2==1 ,sort lpattern(shortdash) lcolor(gs1)) ///
	(line y3 xvar if tag_y3==1 , sort lcolor(gs1) lpattern(dash)), ///
	xlabel(#5) legend(on order(1 "`name1'" 2 "`name2'" 3 "`name3'" )) scheme(s1mono) ///
  	ylabel(#5, angle(horizontal) grid)
	graph save utility_hara_line, replace
***************drop all generated variables*********************************
restore
end


var3_shm  firm_fe_5 firm_fe_w5 firm_fe_m5 wp_shm 20 "wp_size>=5" "All wages" "Women's Wages'" "Men's Wages'" 
graph save Figure_A12, replace

*Generate Figure A19
binscatter firm_fe_5 firm_fe_1 firm_fe_15 workpl_m if main_sample==1, nbins(100) line(none)
graph save Figure_A19L, replace
binscatter firm_fei_5 firm_fei_1 firm_fei_15 workpl_m if main_sample==1, nbins(100) line(none)
graph save Figure_A19R, replace

*Figure A20
 xi: reg educ_year ln_wage flextid if  main==1
predict edu_flex  if  main==1
reg educ_year ln_wage hara_comb if main==1
predict edu_hara if  main==1
 
xi: reg educ_year ln_wage flextid if woman==0  & main==1
predict edu_flextid_m  if woman==0  & main==1
xi: reg educ_year ln_wage flextid if woman==1 & main==1
predict edu_flextid_w  if woman==1& main==1
egen edu_flextid=rowmax( edu_flextid_m edu_flextid_w) if  main==1
reg educ_year ln_wage hara_comb if woman==0 & main==1 & age_cat==3
predict edu_hara_m if woman==0 & main==1 & age_cat==3
xi: reg educ_year ln_wage hara_comb if woman==1 & age_cat==3
predict edu_hara_w if woman==1 & main==1 & age_cat==3
egen edu_hara_comb=rowmax( edu_hara_m edu_hara_w) if age_cat==3

xi: reg educ_year ln_wage part_time if woman==0  & main==1
predict edu_part_m  if woman==0  & main==1
xi: reg educ_year ln_wage part_time if woman==1 & main==1
predict edu_part_w  if woman==1 & main==1
egen edu_part=rowmax( edu_part_m edu_part_w) if  main==1

reg ln_wage `treat' edu_hara_w
reg ln_wage flex edu_flex_w

reg ln_wage `treat' edu_hara_m 
reg ln_wage flex edu_flex_m

bysort  woman: reg ln_wage flex edu_flex
bysort woman: reg ln_wage hara_comb edu_hara_comb


foreach var in  ln_wage {
	foreach sample in  man_main woman_main{

	foreach treat in  hara_comb flextid {
	preserve
	local s=0
	matrix C =J(1,3,.)
	local s=`s'+1
	reg `var' `treat'  [pw=gen_edu2_alder5_w ] if `sample' ==1
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])
	svmat C
	ren C1 est_`s'
	ren C2 ci_low_`s'
	ren C3 ci_high_`s'
	gen pos_`s'=7 if est_`s'!=.
	
	*****
	matrix C =J(1,3,.)
	local s=`s'+1
	xi: reg `var' `treat' i.age_cat i.edu3 i.famtyp europe non_europe    [pw=gen_edu2_alder5_w ] if `sample' ==1 , 
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'
	ren C2 ci_low_`s'
	ren C3 ci_high_`s'
	gen pos_`s'=6 if est_`s'!=.
	
	****
	local s=`s'+1
	xi: reg `var' `treat' firm_fe_5 [pw=gen_edu2_alder5_w ] if `sample' ==1 ,
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'
	ren C2 ci_low_`s'
	ren C3 ci_high_`s'
	gen pos_`s'=5 if est_`s'!=.
	
	local s=`s'+1
	xi: reghdfe `var' `treat'  [pw=gen_edu2_alder5_w ] if `sample' ==1 , abs( sektor ind_code)
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'
	ren C2 ci_low_`s'
	ren C3 ci_high_`s'
	gen pos_`s'=4 if est_`s'!=.
	****
	local s=`s'+1
	xi: reghdfe `var' `treat'    [pw=gen_edu2_alder5_w ] if `sample' ==1 , abs(ssyk3 )
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'
	ren C2 ci_low_`s'
	ren C3 ci_high_`s'
	gen pos_`s'=3 if est_`s'!=.
	****
	local s=`s'+1
	xi: reghdfe `var' `treat' i.age_cat i.edu3 i.famtyp europe non_europe  firm_fe_5    [pw=gen_edu2_alder5_w ] if `sample' ==1 , abs(ssyk3 ind_code sektor)
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'
	ren C2 ci_low_`s'
	ren C3 ci_high_`s'
	gen pos_`s'=2 if est_`s'!=.

	****
	local s=`s'+1
	xi: reg `var' `treat' edu_`treat'  if `sample' ==1  [pw=gen_edu2_alder5_w ] , 
	matrix C[1,1]=_b[ `treat']
	matrix C[1,2]=(_b[ `treat']+1.96*_se[ `treat'])
	matrix C[1,3]=(_b[ `treat']-1.96*_se[ `treat'])

	svmat C
	ren C1 est_`s'
	ren C2 ci_low_`s'
	ren C3 ci_high_`s'
	gen pos_`s'=0 if est_`s'!=.

twoway (scatter  pos_1 est_1  , msize(medlarge)  msymbol(circle)   mcolor(gs4)) ///
	(rbar ci_low_1 ci_high_1 pos_1 , lwidth(vvthin)  fcolor(gs4) lcolor(gs4) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_2 est_2  , msize(medlarge)  msymbol(circle)  mcolor(gs6)) ///
	(rbar ci_low_2 ci_high_2 pos_2 , lwidth(vvthin)  fcolor(gs6) lcolor(gs6) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_3 est_3  , msize(medlarge)  msymbol(circle)  mcolor(gs8)) ///
	(rbar ci_low_3 ci_high_3 pos_3 , lwidth(vvthin)  fcolor(gs8) lcolor(gs8) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_4 est_4  , msize(medlarge)  msymbol(circle)  mcolor(gs10)) ///
	(rbar ci_low_4 ci_high_4 pos_4 , lwidth(vvthin)  fcolor(gs10) lcolor(gs10) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_5 est_5  , msize(medlarge)  msymbol(circle)  mcolor(gs12)) ///
	(rbar ci_low_5 ci_high_5 pos_5 , lwidth(vvthin)  fcolor(gs12) lcolor(gs12) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_6 est_6  , msize(medlarge)  msymbol(circle)  mcolor(gs11)) ///
	(rbar ci_low_6 ci_high_6 pos_6 , lwidth(vvthin)  fcolor(gs14) lcolor(gs14) lpattern(solid) horizontal barwidth(.015)) ///
	(scatter  pos_7 est_7  , msize(medlarge)  msymbol(circle)  mcolor(gs1)) ///
	(rbar ci_low_7 ci_high_7 pos_7 , lwidth(vvthin)  fcolor(gs1) lcolor(gs1) lpattern(solid) horizontal barwidth(.015)) ///
	,xline(0, lcolor(gs7) lpattern(dash)) yscale(range(-1 7))  ///
	 ytitle(, size(zero))  ylabel(none)  xlabel(, grid) scale(.7) ysize(20) xsize(20) scheme(s1mono) title(`sample' `var') ///
	 legend(on order(1 "No controls" 3 "Demographic Controls" 5 "Workplace Pay" ///
	 7 "Sector and Industry Fixed Effects" 9 "Occupation Fixed Effects" 11 "All Controls" ///
	 13 "Bell Compensating Pay") col(4)) 
	graph save reg_`sample'_`var'_`treat', replace
	restore
	}
	}
}
graph combine reg_woman_main_ln_wage_hara_comb.gph  reg_man_main_ln_wage_hara_comb.gph ///
	reg_woman_main_ln_wage_flextid.gph  reg_man_main_ln_wage_flextid.gph, ///
	ysize(18) xsize(20) col(2) row(2) iscale(1) scheme(s1mono) title(Women Men)
graph save Figure A20, replace
		

foreach var in    wpch_ee     { 
	foreach treat in   hara_comb  { 
	foreach samp in man woman{
	preserve
	keep if age <60
	matrix C =J(17,16,.)
	forvalues n=1(1)3 {
		xi: areg `var'_dn`n' `treat'    if `samp'==1 [pw=gen_edu2_alder5_w ], abs(year)
		scalar year= -`n'-.15
		matrix C[`n',1]=_b[`treat']
		matrix C[`n',2]=(_b[`treat']+1.96*_se[`treat'])
		matrix C[`n',3]=(_b[`treat']-1.96*_se[`treat'])
		matrix C[`n',4]=year
		xi: areg `var'_dn`n' `treat'    i.age_cat i.edu3 i.famtyp europe non_europe      if `samp'==1 [pw=gen_edu2_alder5_w ], abs(year)
		scalar year= -`n'
		matrix C[`n',5]=_b[`treat']
		matrix C[`n',6]=(_b[`treat']+1.96*_se[`treat'])
		matrix C[`n',7]=(_b[`treat']-1.96*_se[`treat'])
		matrix C[`n',8]=year
		xi: areg `var'_dn`n' `treat'  i.age_cat  i.edu3 i.famtyp europe non_europe workpl_m   if `samp'==1 [pw=gen_edu2_alder5_w ], abs(year)
		scalar year= -`n'+.15
		matrix C[`n',9]=_b[`treat']
		matrix C[`n',10]=(_b[`treat']+1.96*_se[`treat'])
		matrix C[`n',11]=(_b[`treat']-1.96*_se[`treat'])
		matrix C[`n',12]=year
	}
	forvalues n=1(1)5 {
		local v= `n'+4		
		xi: areg `var'_d`n' `treat'    if `samp'==1, abs(year)
		matrix C[`v',1]=_b[`treat']
		matrix C[`v',2]=(_b[`treat']+1.96*_se[`treat'])
		matrix C[`v',3]=(_b[`treat']-1.96*_se[`treat'])
		matrix C[`v',4]=year
		xi: areg `var'_d`n' `treat'   i.age_cat i.edu3 i.famtyp europe non_europe    if `samp'==1 [pw=gen_edu2_alder5_w ],  abs(year)
		scalar year= `n'
		matrix C[`v',5]=_b[`treat']
		matrix C[`v',6]=(_b[`treat']+1.96*_se[`treat'])
		matrix C[`v',7]=(_b[`treat']-1.96*_se[`treat'])
		matrix C[`v',8]=year
		xi: areg `var'_d`n' `treat'  i.age_cat fem_bos i.edu3 i.famtyp europe non_europe workpl_m  if `samp'==1 [pw=gen_edu2_alder5_w ], abs(year)
		scalar year= `n'+.15
		matrix C[`v',9]=_b[`treat']
		matrix C[`v',10]=(_b[`treat']+1.96*_se[`treat'])
		matrix C[`v',11]=(_b[`treat']-1.96*_se[`treat'])
		matrix C[`v',12]=year
	}
	scalar year= -.2
	matrix C[4,1]=0
	matrix C[4,4]=year
	scalar year= 0
	matrix C[4,5]=0
	matrix C[4,8]=year
	scalar year= 0.2
	matrix C[4,9]=0
	matrix C[4,12]=year
	svmat C
	rename C1 a_estimate
	renae C2 a_upper
	rename C3 a_lower
	rename C4 a_year
	rename C5 o_estimate
	rename C6 o_upper
	rename C7 o_lower
	rename C8 o_year
	rename C9 e_estimate
	rename C10 e_upper
	rename C11 e_lower
	rename C12 e_year
	
	twoway (rbar a_upper a_lower a_year if a_upper!=a_lower, fcolor(black) lcolor(none) barwidth(.05)) ///
		(scatter a_estimate a_year , mcolor(black) msize(medium) msymbol(circle)) ///
		(rbar o_upper o_lower o_year if a_upper!=o_lower , fcolor(gs4) lcolor(none) barwidth(.05)) ///
		(scatter o_estimate o_year , mcolor(gs4) msize(medium) msymbol(square)) ///
		(rbar e_upper e_lower e_year if e_upper!=e_lower , fcolor(gs8) lcolor(none) barwidth(.05)) ///
		(scatter e_estimate e_year , mcolor(gs8) msize(medium) msymbol(triangle)) ///
        , yline(0, lwidth(medthick) lpattern(dash) lcolor(black))  xline(-.5 .5, lcolor(black) lpattern(dot) lwidth(medthick)) ///
		title(`samp', size(medium)) xtitle(Years from measurement)  xlabel(-3(1)5) ylabel(  , angle(horizontal)grid) ///
		ytitle(Estimated difference SH-no SH) scheme(s1mono) ///
		legend(order(2 "Survey Year Fixed Effects" 4 " Demographic Controls" ///
		6 "Workplace Share of Men") col(1) region(lcolor(none)))
		drop a_estimate-e_year
		graph save event_`var'_`samp', replace
		restore
	
	}
}
}
graph combine event_wpch_ee_woman.gph   event_wpch_ee_man.gph  , col(2) scheme(s1mono)  ysize(7) xsize(20)  iscale(1.4)
graph save Figure_6, replace
